Select between five and 30 locations within a city or metropolitan area. For each location, generate isochrones for the same travel time for at least two different modes. Calculate the area of each isochrone and compare the areas of the isochrones for the two (or more) modes you analyzed. Create three figures (which may or may not be maps) to illustrate the results of your analysis.
library(osmdata)
library(opentripplanner)
library(tidyverse)
library(sf)
library(ggthemes)
library(ggspatial)
library(ggmap)
library(stringr)
library(sp)
library(rgeos)
library(tidygeocoder)
library(raster)
childcareaddress <- read.csv("2miradiuschildcare.csv")
childcareaddress <- as.data.frame(childcareaddress)
opq(bbox = 'Boston MA USA') %>%
add_osm_feature(key = 'highway') %>%
osmdata_xml(file = 'OTP/graphs/default/boston_streets.osm')
MA_state_plane <- "+proj=lcc +lat_1=41.71666666666667 +lat_2=42.68333333333333 +lat_0=41 +lon_0=-71.5 +x_0=200000 +y_0=750000 +ellps=GRS80 +units=m +no_defs"
boston_street_features <- opq(bbox = 'Boston MA USA') %>%
add_osm_feature(key = 'highway') %>%
osmdata_sf()
boston_streets <- boston_street_features$osm_lines %>%
st_transform(crs = MA_state_plane)
path_otp <- otp_dl_jar("OTP")
otp_setup(otp = path_otp, dir = path_data, memory =1024)
otpcon <- otp_connect()
## Router http://localhost:8080/otp/routers/default exists
I found this list of addresses of certified childcare centers in Mattapan from the Mass.gov website. I cleaned up the data in excel to comport with geocoding the addresses using tidygeocoder.
address_list = c("35 Astoria Street, Mattapan, MA",
"40 HOSMER STREET, Mattapan, MA",
"32 WOODBOLE AVENUE, Mattapan, MA",
"7 LANDOR ROAD # 1, Mattapan, MA",
"23 WILMORE STREET #1, Mattapan, MA",
"5 Mildred Ave, Mattapan, MA",
"500 Walk Hill Street, Mattapan, MA",
"97 WEST SELDEN ST, Mattapan, MA",
"608 NORFOLK STREET, Mattapan, MA",
"130 Orlando Street, Mattapan, MA",
"62 RIDGEVIEW AVE, Mattapan, MA",
"55 Hollowell Street, Mattapan, MA",
"760 Morton Street #5, Mattapan, MA",
"100 ROCKDALE ST, Mattapan, MA",
"10 ELLISON AVE, Mattapan, MA",
"6 Lorna Rd, Mattapan, MA",
"1295 BLUE HILL AVE, Mattapan, MA",
"14 Groveland Street #1, Mattapan, MA",
"70 FAVRE ST # 2, Mattapan, MA",
"124 LORNA ROAD, Mattapan, MA",
"37 Duke Street #2, Mattapan, MA",
"14 DANIA STREET, Mattapan, MA",
"37 Babson Street, Mattapan, MA",
"5 Mildred Ave, Mattapan, MA",
"778 MORTON STREET # 1, Mattapan, MA",
"624 HARVARD STREET, Mattapan, MA",
"438 River Street, Mattapan, MA",
"535 RIVER STREET, Mattapan, MA",
"8 LESTON ST, Mattapan, MA",
"596 RIVER ST, Mattapan, MA",
"5 MILDRED AVENUE, Mattapan, MA",
"22 HARMON STREET, Mattapan, MA",
"73 WOODBOLE AVENUE #216, Mattapan, MA",
"74 WESTMORE ROAD, Mattapan, MA",
"5 Mildred Avenue, Mattapan, MA",
"712 MORTON STREET #3, Mattapan, MA",
"119 Lorna Road, Mattapan, MA",
"36 Rockingham Road, Mattapan, MA",
"29 Woodruff Way, Mattapan, MA",
"9 GREENDALE RD # 2, Mattapan, MA",
"54 Hiawataha Road, Mattapan, MA",
"1009 Morton Street, Mattapan, MA",
"163 COLORADO ST, Mattapan, MA",
"255 RIVER ST, Mattapan, MA",
"152 STANDARD STREET # 129, Mattapan, MA",
"55 Woodhaven Street Apt. 2, Mattapan, MA",
"12 ELLISON AVENUE # 1, Mattapan, MA",
"557 NORFOLK ST, Mattapan, MA",
"100 Hebron Street, Mattapan, MA",
"20 Outlook Road, Mattapan, MA")
points <- geo(address = address_list, mode = "batch")
head(points)
## # A tibble: 6 x 3
## address lat long
## <chr> <dbl> <dbl>
## 1 35 Astoria Street, Mattapan, MA 42.3 -71.1
## 2 40 HOSMER STREET, Mattapan, MA 42.3 -71.1
## 3 32 WOODBOLE AVENUE, Mattapan, MA 42.3 -71.1
## 4 7 LANDOR ROAD # 1, Mattapan, MA 42.3 -71.1
## 5 23 WILMORE STREET #1, Mattapan, MA 42.3 -71.1
## 6 5 Mildred Ave, Mattapan, MA 42.3 -71.1
points2 <- na.omit(points)
points2 <- st_as_sf(x = points2,
coords = c("long", "lat"),
crs = 4326)
I then re-merged the geocoded addreses with the original csv in order to not lose the data associated with the addresses.
merged.data <- merge(points2, childcareaddress, by="address")
With my geocoded childcare addresses ready to go, I calculated the 5 minute walking and biking isochrones around each address.
multiple_points_5min_walk <-
otp_isochrone(otpcon = otpcon, fromPlace = points2,
mode = "WALK", cutoffSec = 300)
multiple_points_5min_bike <-
otp_isochrone(otpcon = otpcon, fromPlace = points2,
mode = "BICYCLE", cutoffSec = 300)
ggplot(multiple_points_5min_walk) +
annotation_map_tile(zoomin = 1, progress = "none") +
geom_sf(fill ="blue", alpha=0.2) +
theme_map()
## although coordinates are longitude/latitude, st_intersects assumes that they are planar
## although coordinates are longitude/latitude, st_intersects assumes that they are planar
bus_stops <- st_read("https://opendata.arcgis.com/datasets/55586c8f54954f8e8fae5f40cb953d15_0.kml?outSR=%7B%22latestWkid%22%3A26986%2C%22wkid%22%3A26986%7D",
quiet=TRUE)
multiple_points_5min_walk <- multiple_points_5min_walk %>%
mutate(transit_score = lengths(st_covers(geometry, bus_stops)))
ggplot(multiple_points_5min_walk) +
annotation_map_tile(zoomin = 1, progress = "none") +
geom_sf(aes(fill=transit_score), alpha=.5) +
theme_map()
## although coordinates are longitude/latitude, st_intersects assumes that they are planar
## although coordinates are longitude/latitude, st_intersects assumes that they are planar
iso_5min_walk <-
otp_isochrone(otpcon = otpcon, fromPlace = points2,
mode = "WALK", cutoffSec = 300) %>%
st_transform(crs = MA_state_plane) %>%
mutate(mode = "walk")
iso_5min_bike <-
otp_isochrone(otpcon = otpcon, fromPlace = points2,
mode = "BICYCLE", cutoffSec = 300) %>%
st_transform(crs = MA_state_plane) %>%
mutate(mode = "bicycle")
iso_all_modes <- rbind(iso_5min_walk, iso_5min_bike)
ggplot(iso_all_modes) +
annotation_map_tile(zoomin = 2, type = "cartolight", progress = "none") +
geom_sf(aes(fill = mode), alpha = 0.5) +
geom_sf(data = points2) +
scale_fill_viridis_d(name = "Area that is reachable within 5 minutes",
labels = c("By bike", "By foot")) +
theme_map() +
labs(caption = "Basemap Copyright OpenStreetMap contributors")
iso_areas <- iso_all_modes %>%
mutate(area = st_area(iso_all_modes)) %>%
st_set_geometry(NULL) %>%
pivot_wider(names_from = mode, values_from = area)
iso_areas <- iso_areas %>%
filter(bicycle != "NULL", walk != "NULL") %>%
filter(str_detect(bicycle,"c")==FALSE)
ggplot(iso_areas,
aes(x = as.numeric(walk), y = as.numeric(bicycle))) +
geom_point() +
scale_x_continuous(name =
"Area within a five-minute walking distance\nof a childcare center\n(square km)",
breaks = breaks <- seq(10000, 260000, by = 20000),
labels = breaks / 1000000) +
scale_y_continuous(name =
"Area within a five-minute driving distance\nof a childcare center\n(square km)",
breaks = breaks <- seq(0, 2000000, by = 100000),
labels = breaks / 1000000) +
theme_bw()
mergedgeometry <- read.csv("mergeddatafromPlace.csv")
mergecapacityiso <- merge(mergedgeometry, iso_areas, by="fromPlace")
ggplot(mergecapacityiso,
aes(x = as.numeric(walk), y = as.numeric(bicycle), size = Capacity)) +
geom_point() +
scale_x_continuous(name =
"Area within a five-minute walking distance\nof a childcare center\n(square km)",
breaks = breaks <- seq(10000, 260000, by = 20000),
labels = breaks / 1000000) +
scale_y_continuous(name =
"Area within a five-minute biking distance\nof a childcare center\n(square km)",
breaks = breaks <- seq(0, 2000000, by = 100000),
labels = breaks / 1000000) +
scale_size_continuous(name = "Childcare Center Capacity")+
theme_bw()